###########################################
# Table 1: Summary statistics

###########################################
# Select variables to present in summary statistics
sumstat <- df[,c("Pass","perc_yes","itm_prop_change","post","presidential_el","midterm_el","odd_year","turnout","failed_recently","cost_bond_ph","levy_increase")]

# 
sumstat$itm_prop_change <- sumstat$itm_prop_change*100
sumstat$turnout <- sumstat$turnout*100

# Rename for readibility
colnames(sumstat) <- c("Passed","Approval Rates (%)","ChangeDed (%)","Post","Presidential election","Midterm election","Odd year election","Voters' turnout (%)","Recently defeated","Cost bond per housing unit ( )","Parcel levy increase (%)")

# Separate between pre and post TCJA
sumstat_pre <- sumstat[sumstat$Post ==0,]
sumstat_post <- sumstat[sumstat$Post ==1,]
sumstat$Post <- sumstat_post$Post <- sumstat_pre$Post <- NULL


##
# Compute selected summary statistics
##

results <- matrix(nrow = ncol(sumstat),ncol = 8)
for(i in 1:ncol(sumstat_post)){
  results[i,1] <- mean(sumstat[,i],na.rm = TRUE) # mean 
  results[i,2] <- sd(sumstat[,i],na.rm = TRUE) # standard deviation
  results[i,3] <- median(sumstat[,i],na.rm = TRUE) # median
  
  results[i,4] <- mean(sumstat_pre[,i],na.rm = TRUE) # mean pre-TCJA
  results[i,5] <- mean(sumstat_post[,i],na.rm = TRUE) # mean post-TCJA
  a <- t.test(sumstat_pre[,i], sumstat_post[,i], alternative = "two.sided", var.equal = FALSE) # Test difference in means
  results[i,6] <- results[i,4] - results[i,5] # difference in means
  results[i,7] <- a$statistic # t-statistics of difference in means
  results[i,8] <- a$p.value # p-values of difference in means
  
}

##
# Create the output
##

results <- as.data.frame(results)
rownames(results) <- colnames(sumstat)
colnames(results) <- c("Mean","Std. dev.","Median","Mean - Pre","Mean - Post","Difference","t-statistics","pvalues")
results$stars <- ifelse(results$pvalues<0.01,3,ifelse(results$pvalues<0.05,2,ifelse(results$pvalues<0.1,1,0)))

stargazer(results[,c(1:7,9)], 
          summary = FALSE,
          type = "latex",
          digits = 2,
          out = "Output/Table1.Summary Statistics.tex",
          title = "Summary statistics",
          label = "tab:refdum_sumstat")


# Clean workspace
rm(list = setdiff(ls(),c("df","city")))
